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ABSTRACT 


A review is made of the performance of a variety of turbulence models in the 
evaluation of a particular well-documented transonic flow. This is done to supple- 
ment a previous attempt to calibrate and verify transonic airfoil codes by including 
many more turbulence models than used in the earlier work and applying the calcula- 
tions to an experiment that did not suffer from uncertainties in angle of attack and 
was free of wind tunnel interference. It is found from this work, as well as in the 
earlier study, that the Johnson-King turbulence model is superior for transonic 
flows over simple aerodynamic surfaces, including moderate separation. It is also 
shown that some field equation models with wall function boundary conditions can be 
competitive with it. 


INTRODUCTION 


The accuracy of aerodynamic computations, at realistic Reynolds numbers where 
boundary-layer turbulence exists, depends on two distinct, yet interrelated ele- 
ments. These are the numerical and the physical aspects of the solution scheme. 
Currently, practical aerodynamic problems have to be solved in a manner where the 
physical aspects of turbulence, present at the Reynolds numbers of interest, are 
treated with statistical models. Although such turbulence models can be avoided 
with time- and space-accurate solutions of the dynamics of turbulence, to date such 
computations have been confined to very simple flow fields at extremely low Reynolds 
numbers. Their costs are much too high to be considered for design work in the 
foreseeable future (ref.1), even if sub-grid models can be developed to allow han- 
dling the Reynolds numbers of practical interest. 

The statistical equations for turbulent flows, corresponding to the basic 
Navier-Stokes equations for compressible flows, are derived by first expanding the 
instantaneous values of the dependent variables, such as velocity components, den- 
sity, and temperature into the sum of mean and fluctuating quantities. The equa- 
tions are then either averaged in time, over periods much longer than the time 
scales of the largest turbulence eddies, or are ensemble-averaged. When transonic 
flows are considered, simpler equations result when the velocity terms are weighted 
with the local density in the averaging process. Except for existence or absence of 
the time-dependent terms, the resulting equations have the same form for either 
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time- or ensemble-averaging and are called the Reynolds-averaged equations. The 
dependent variables in these equations are composed of the mean quantities of the 
velocity and state variables and of moments of their instantaneous quantities. 

These moments are called Reynolds stresses or heat flux, and it is the evaluation of 
these quantities that is central to physical turbulence modeling. 

In recent years, there has been a considerable expenditure of man-power and 
computer costs in developing numerical techniques for solving the Reynolds-averaged 
equations. Today, the flow fields about three-dimensional bodies as complicated as 
wing-body configurations can be computed in a reasonable time with numerically 
accurate techniques. The basic algorithms are becoming more efficient and robust, 
complex meshes are being generated routinely, and codes run sufficiently economi- 
cally so that mesh-independent solutions can be tested and demonstrated. Color 
graphics permit detailed and illuminating visualization of the results. Yet despite 
this enormous rate of progress and maturation of the field of computational fluid 
dynamics (CFD), the design community still regards the results to be useful only in 
a qualitative sense and is pressuring the CFD community to calibrate and verify 
their codes against carefully obtained experimental data so that these codes can be 
established as accurate design tools. 

Because the numerical aspects of a problem can be systematically isolated and 
improved upon, differences between the output of a good and carefully applied numer- 
ical code and experimental data are generally identified as a weakness of the turbu- 
lence model. Of course, not all computations are done with the proper care, so that 
some differences between experimental data and computed results can be numerical as 
well as physical. That this occurs in published results will be demonstrated by 
some of the examples shown in the following section. 

One attempt to calibrate and verify some existing codes for computing two- 
dimensional airfoils at transonic speeds was the Viscous Transonic Airfoil Workshop 
held at the AIAA 25th Aerospace Sciences Meeting at Reno, Nevada, in January 1987, 
and reported in reference 2. Fifteen different author groups participated, with 
codes that included 16 Navier-Stokes methods, 2 Euler/boundary-layer methods, and 5 
full-potential/boundary-layer methods. The turbulence models variously employed in 
the Navier-Stokes methods included the classical algebraic eddy viscosity model of 
Cebeci-Smith (ref. 3), the Baldwin-Lomax algebraic eddy viscosity model (ref. 4), 
the Johnson-King model (ref. 5), which accounts for turbulence lag with an ordinary 
differential equation for growth of the maximum shear stress that is then used to 
scale an algebraic eddy viscosity , and the Coakley q-omega model (ref. 6), the sole 
representative of the class of two field equation models that define length and 
velocity scales for an eddy-viscosity evaluation throughout the boundary layer, 
including the region near the surface. The interactive inviscid-/viscous-layer 
codes either use the Cebeci-Smith model in the differential equations applying to 
the boundary layer or an integral model. 

The results of this workshop will be summarized in the following section. 

These results have to be considered to be tentative because they apply to airfoil 
experiments that experienced wall effects, with an attendant uncertainty in the 
experimental data, and, also, only a limited number of turbulence models were 
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represented in the Navier-Stokes computations. To complement this work, a later 
section of this paper is devoted to the detailed comparison of computed results, 
from a single computer code that contains a larger number of turbulence models, with 
the data of an experiment that is relatively free of wind tunnel wall effects and 
contains measurements of the turbulence itself. These comparisons, as well as the 
earlier work with the airfoils, are then used here to develop some conclusions 
regarding the state of the art of modeling turbulence in two-dimensional transonic 
flows. 


SUMMARY OF VISCOUS TRANSONIC AIRFOIL WORKSHOP 


The results of the AIAA Viscous Transonic Airfoil Workshop (ref. 2) are summa- 
rized on figure 1, where the coefficient of lift is plotted against the coefficient 
of drag, expressed in counts, for three different airfoils at various Mach numbers 
and angles of attack. On this figure, the solid circles represent the experimental 
data. The open circles represent the inviscid/viscous interaction schemes, without 
distinction as to whether Euler or full potential equations apply to the inviscid 
region or whether differential or integral methods are used in the boundary layer. 
The remaining symbols refer to Navier-Stokes calculations. The open squares repre- 
sent the Baldwin-Lomax algebraic eddy viscosity model; the squares with ticks repre- 
sent the Cebeci-Smith model; the diamonds represent the Johnson-King lag model, and 
finally, the triangle represents Coakley's two-equation model. 

Figure 1(a) shows the comparison of calculated CL(CD) with the experimental 
data on an NACA 0012 airfoil (ref. 7) for three conditions of Mach number and angle 
of attack, all at a chord Reynolds number of nine million. In the lower left-hand 
corner, there are shown the data and computed results at M = 0.7 and at an angle 
of attack of 1.49°, which includes wind tunnel interference corrections recommended 
by the experimenter. All of the computed results yield results that are within 10% 
of the lift, and only two of the inviscid/viscous interactive methods miss the drag 
by as much as 20 counts. Thus, for attached flows that are weakly transonic, all 
the computation methods give reasonable quantitative results. 

The datum point and computed values of Cp(Cp), at M = 0.55 and a corrected 
angle of attack of 8.34°, are shown in the upper portion of figure 1(a). Under 
these conditions, the flow over the top of the airfoil had a supersonic bubble 
followed by a shock wave at about 0.1 chord, where a slight separation occurred. 
Again, except for two of the viscous/inviscid interactive schemes, the lift is 
calculated to be between 1% high and Q% low by the various methods. The closest 
points to the experimental values in terms of Cp and Cp are the Navier-Stokes 
solutions using the Baldwin— Lomax and Cebeci-Smith models. It is important to note 


583 


that equivalent methods, i.e. those solving, ostensibly, the same conservation and 
modeling equations, also show significant differences, especially in drag. For 
example, the two Navier-Stokes codes that utilized the Johnson-King model, the 
diamonds, show differences of over 50 counts in drag coefficient, although lift only 
differs by about Similar differences are seen in the results produced by codes 

employing the Baldwin-Lomax or Cebeci-Smith models. In this example, it is seen 
that the numerical differences between codes with supposedly the same turbulence 
models can be larger than the differences produced by the different turbulence 
models within a single code. These observations bring into serious question the 
premise stated earlier that code verification processes will be tests primarily of 
the physical, rather than the numerical aspects of a code. It appears, in this 
example, that the numerical aspects of some of the codes may have not been treated 
with sufficient care to make the solutions of each code numerically "correct”, 
thereby permitting their solutions to be differentiated by their physical modeling. 

A point that has to be made here is that there are often some arbitrary deci- 
sions made in the application of a particular turbulence model, and that some of the 
differences shown in the last example could be caused by differences in how those 
decisions were made rather than related to the numerics. For example, in the 
Cebeci-Smith model, the displacement thickness is used to define the eddy viscosity 
in the outer region of the boundary layer. The evaluation of a displacement thick- 
ness is unambiguous only if the irrotational flow outside the boundary layer pro- 
duces a uniform velocity component parallel to the surface. This is not the case in 
the current example and in most complex flows. In fact, this difficulty in defining 
the displacement thickness was the underlying reason for the changes that were 
introduced into the Baldwin-Lomax model. This latter model also can suffer from 
ambiguities, as its velocity scale in the outer regions of the boundary layer can be 
multivalued. This results in an imprecise prescription of the eddy viscosity in 
this region. Another significant choice that can introduce differences in the 
solutions can be the distance of the first mesh point from the surface. If this 
point does not lie well within the viscous-dominated region of the boundary layer, 
errors in the skin friction determined by a particular turbulence model can occur. 
From the published results, it is not possible to discern what differences such 
choices may have made in the various solutions shown on the figure. 

A third example of flow over a NACA 0012 airfoil is shown in the lower right of 
figure 1(a). Here the Mach number is 0.799 and the corrected angle of attack is 
2.26°. Under these conditions, a shock wave developed on the upper surface at about 
mid-chord that was sufficiently strong to cause significant boundary-layer separa- 
tion. Again, ostensibly similar methods yield a spread in the results. In this 
case, however, the two computations with the Johnson-King model agree rather well 
and, in terms of both and Cq, are the closest to the experimental result, except 
for the same interactive method that gave the best results in the previous example 
(ref. 8). In this case the two-equation model gave poorer results than some of the 
algebraic eddy viscosity models. Also, as in the previous example, most of the 
codes gave values of drag much larger than the experimental value. 
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Figure 1(b) shows C L (C D ) for the RAE 2822 airfoil for two test conditions in 
reference 9. The symbols have the same meaning as in figure 1(a). Both sets of 
experimental conditions were obtained at chord Reynolds numbers of 6.5*10^. The 
data point on the left represents the test results at M = 0.725 at a corrected 
angle of attack equal to 2.92°. A shock wave developed on the upper surface at 
about X/C = 0.52; but its strength was insufficient to separate the boundary 
layer. Under these conditions, the bulk of the computational results cluster about 
the experimental datum point within ±10$ on either axis. It is surprising that the 
computation with the Johnson-King model, the diamond, which was generally the more 
successful model in figure 1(a), yields one of the poorer results in this case. 

Here too, the differences between supposedly similar methods show differences as 
large as those identified with different turbulence models, albeit both sets of 
differences are relatively small. 

The point on the right of this figure corresponds to test conditions at higher 
Mach number and angle of attack, namely, M = 0.75 and 3.19°. Under these wind 
tunnel conditions, the upper surface experiences a shock wave, still at X/C = 0.52, 
but with a strength sufficient to cause boundary-layer separation according to most 
of the turbulence models. Here the better computed results are due to three 
viscous/inviscid interactive models, including that of reference 8, and the 
Johnson-King model in a Navier-Stokes code. The Coakley two-equation model and the 
eddy viscosity models cluster together at high values of drag. 

The third figure of this group, figure 1(c), relates to an experiment (ref. 10) 
conducted with a supercritical transonic flow airfoil, designated the VA-2, that was 
developed in West Germany. For this case, the only computed results that exist are 
from a single code developed by Coakley (ref. 10) which contains the Cebeci-Smith, 
Johnson-king , and Coakley' s q-omega models. In this case, the tests were conducted 
at a chord Reynolds number of 6*1 0^ and a wind tunnel angle of attack of 1°. The 
point on the left of the figure represents the datum at M = 0.73, where a shock 
wave developed at X/C = 0.35. All the models overpredict the lift, and underpre- 
dict the drag. As with most of the cases with the other airfoils, the two-equation 
model gave results essentially identical with those of the Cebeci-Smith model. The 
Johnson-King model yielded lift within a few percent of the experimental value, but 
the drag was 30 counts too low, and no better than the values from the other models. 

The point on the right of figure 1(c). corresponds to a Mach number of 0.78, 
where the shock wave moved to X/C = 0.645, rather far back on the airfoil. Again, 
the lift is predicted well by the Johnson-King model, whereas the other models do 
much better with regard to drag. 

Several conclusions can be drawn from the above discussion of the comparison of 
calculated C^(Cp) with the experimental data from three different airfoils in 
transonic flow. These are as follows: 

1. In two-dimensions, some interactive inviscid/viscous methods can be competi- 
tive with Navier-Stokes computations containing the better turbulence models. 
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2. In the Navier-Stokes codes, the algebraic models of Cebeci-Smith and 
Baldwin-Lomax yield similar results, which are generally equivalent to results from 
the the more complex Coakley two-equation q-omega model. 

3. The Navier-Stokes codes that yield the best results in most of the cases 
employ the Johnson-King model. 

These must be treated as tentative conclusions. Although the data represent 
the state of the art, they suffer, as do most wind tunnel experiments, from wall 
interference effects, especially those cases where the airfoils experienced larger 
regions of separation. A manner of accounting for the wind tunnel effects is to 
"correct" the geometrical angle of attack to some effective value. In the computa- 
tions cited above, the effective angles of attack used by the various author groups 
differed from each other and the reasons for selecting particular values were not 
stated in reference 2. How much the uncertainty in the angle of attack in each case 
affects the results shown in reference 2 is not clear. In addition, the assessment 
in reference 2 suffers from a limited number of turbulence models, and in the case 
of the two-equation model does not represent the best of that class. 

To provide a broader assessment of the state of two-dimensional transonic 
turbulence modeling, we decided to complement the work described above by computing 
the flow field of a transonic flow whose data have been shown to be free of wind 
tunnel interference (ref. 11) and angle of attack uncertainty, and to include addi- 
tional turbulence models. Since all but one of these models are included in a 
single Navier-Stokes code, the uncertainties introduced by different codes and 
solution algorithms are also essentially eliminated. These results are described 
and discussed in the following section. 


TRANSONIC FLOW TURBULENCE MODEL COMPARISONS 


Experimental Basis for Comparison 

The experiment chosen as the standard for the turbulence model comparisons 
utilized the transonic flow over a circular arc bump affixed to a hollow circular 
cylinder aligned with the flow direction (ref. 11). The bump acts aerodynamically 
in a manner similar to the top surface of an airfoil. The configuration of this 
wind tunnel model and the flow conditions for which most of the turbulence model 
comparisons were made are shown on figure 2. The thin-walled cylinder was 15.2 cm 
in diameter and extended 61 cm upstream of the bump's leading edge. The bump had a 
chord of 20.32 cm and a thickness of 1.91 cm. At its leading edge, the bump 
departed from its circular arc shape and blended smoothly into the upstream cyl- 
inder. Natural transition occurred on the cylinder, and the resulting turbulent 
boundary-layer thickness at the leading edge of the bump was about 1 cm. This was 
sufficiently thick to allow boundary layer measurements, but not so thick as to 
introduce significant streamwise curvature effects. The measurements included 
surface pressure, mean velocities, and the Reynolds stresses. These latter data, 
being profile measurements, allow examination of the reasons for the success or 
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failures of individual turbulence models. Another decided advantage of using this 
experiment as a standard is the fact that the wind tunnel model was run in both the 
NASA Ames Research Center 2x2 ft and 6x6 ft transonic wind tunnels and that 
the pressure distributions at the same Mach and Reynolds numbers were essentially 
the same (ref. 12). This suggests that the data used to test the turbulence models 
are relatively free of wind tunnel interference. Of course, another advantage of 
this data set is the lack of ambiguity regarding the angle of attack, as was the 
case with the airfoils. A shortcoming of using these data, however, is that they 
were used to develop the Johnson-King model and, therefore, it would be expected 
that the best performance other models could be made to achieve would be to approach 
the behavior of this lag model. Therefore, testing a variety of turbulence models 
with these data can lead to a necessary condition for other turbulence models to 
satisfy, but is not sufficient for ascertaining the relative performance of turbu- 
lence models for flows with separation more massive than occurred within the refer- 
enced experiment. 


Numerical Procedure 

The primary numerical procedure used in computations shown here was originally 
developed by C. C. Horstman (ref. 13) based on MacCormack's 1981 algorithm 
(ref. 14). The code has been modified into a series of codes to accommodate the 
various turbulence models compared in this paper. The MacCormack method used here 
is basically an explicit second-order, predictor-corrector finite-difference method, 
modified by an efficient implicit algorithm. In the calculations, the computational 
domain extended in the flow direction from — 1 40 cm to +90 cm, relative to the onset 
of the bump. In the radial direction the domain extended from the surface of the 
model to a radius of 90 cm. The mesh size was 129 x 45, and allowed variable spac- 
ing in both the axial and radial directions. In the streamwise direction, the mesh 
varied from 0.16 cm, near the shock wave, to 12 cm at the downstream boundary. 

Normal to the surface, an exponentially stretched mesh was used that allowed the 
first mesh point from the surface to lie well within the sublayer, for those cases 
where integrations were performed to the surface. When wall functions were used, 
the near wall mesh centers were generally well into the fully turbulent flow, except 
at the stations close to and within the separation zone. 

The upstream boundary conditions were prescribed by uniform free stream condi- 
tions, with the total pressure and temperature held constant and the static proper- 
ties found from relationships based on the method of characteristics. The down- 
stream conditions set all axial gradients to zero. In all cases, the surface 
boundary conditions were zero velocity and near adiabatic wall temperature. Uniform 
free stream conditions were applied at the far-radial boundary, as justified in 
reference 15. 

As discussed in the next section, the Wilcox-Rubesin- turbulence model (ref. 16) 
required a different numerical algorithm, namely MacCormack's 1976 hybrid method, 
reference 17. 
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Turbulence Models Considered in Comparative Analysis 

In this comparison of turbulence models, the following turbulence models were 
applied: Cebeci and Smith (ref. 3), Johnson and King (ref. 5), Jones and Launder 
(ref. 18), Viegas and Rubesin (refs. 19 and 20), Coakley (refs. 6 and 21), and 
Wilcox and Rubesin (ref. 16). One feature characteristic of all these models is 
that Reynolds stresses are expressed in terms of an eddy viscosity in a Boussinesq 
constitutive relationship. 

The Cebeci-Smith turbulence model divides the boundary layer into an inner and 
outer zone, where the eddy viscosity is expressed algebraically in each layer. In 
the inner layer, the eddy viscosity is expressed with a length scale that is the 
distance from the surface and and with a velocity scale that is the shear velocity 
at the wall. In the outer layer these scales are the displacement thickness, com- 
puted neglecting density variations, and the boundary-layer edge velocity. The 
extent of the inner layer is determined by requiring the minimum of the eddy viscos- 
ities corresponding to the inner and outer layers be used at each position from the 
surface. This model has become the standard for use with attached boundary-layer 
flows. In its current application (ref. 15) the van Driest damping function, used 
to account for viscous effects in the sublayer and buffer layer, was corrected for 
non-zero values of the streamwise pressure gradient. 

Next in mathematical complexity is the turbulence model of Johnson and King. 
This model also divides the boundary layer into two zones where the eddy viscosities 
are expressed algebraically. In the inner layer, however, the velocity scale used 
in this model depends on the maximum shear stress within the boundary layer, not the 
wall shear stress (e.g. expressed as shear velocity) as in the Cebeci and Smith 
model. This maximum shear stress is found from the solution of an ordinary differ- 
ential equation that represents the rate of the development of the maximum shear 
stress with distance along the surface. The equation is expressed conveniently in 
terms of the departure of shear stress and eddy viscosity from their corresponding 
quantities in equilibrium with the mean motion. Compared to the Cebeci and Smith 
model, the Johnson and King model contains two additional modeling coefficients, one 
that scales the rate of development of the maximum shear stress and another that 
scales the effect of turbulence diffusion within that rate equation. The values of 
these coefficients were chosen by comparison with the data of reference 18 , as well 
as with those from the bump experiment (ref. 11). The computed results using this 
model are from reference 5. 

The next turbulence model to be discussed is that of Jones and Launder, refer- 
ence 19- This model is an example of turbulence models that employ two additional 
partial differential field equations to evaluate the velocity and length scales of 
the turbulence, thereby allowing for the rate processes in the development in the 
local state of the scales of turbulence. In this model, the dependent variable in 
the first equation is the kinetic energy of turbulence, whereas that in the second 
equation is the rate of dissipation of the turbulence kinetic energy. In establish- 
ing the eddy viscosity, the kinetic energy of turbulence acts as the velocity scal- 
ing parameter, while the length scale is composed of an algebraic combination of the 
kinetic energy and the rate of dissipation of kinetic energy. This model also 
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possesses the feature of special viscous terms in these turbulence field equations 
that allow integrations to the surface. The model has been extended to account for 
compressibility, necessary in the computation of transonic flows, by use of mass 
weighting the turbulence equations and allowing for non-zero mean velocity diver- 
gence in the constitutive relationship for the eddy viscosity. 

It has been shown that the popular Jones-Launder turbulence model yields poor 
values of skin friction in the presence of a shock wave in transonic flow, and that 
the introduction of wall functions, originally designed to save computation costs, 
surprisingly improved the prediction of skin friction (ref. 20). Accordingly, the 
compressible-flow, wall-function models of Viegas and Rubesin were applied to the 
Jones-Launder model for the test case considered here. These wall-function 
approaches allowed use of fewer mesh points normal to the surface and reduced the 
stiffness of the turbulence equations; both effects lead to considerable saving in 
computation costs over that of the Jones-Launder model itself. 

These wall-function methods depend on establishing the wall shear from a knowl- 
edge of the mean velocity at the center of the first mesh volume off the surface and 
from the turbulence kinetic energy at the centers of the first two points off the 
surface. The latter kinetic energy information reflects the local axial pressure 
gradient along the surface. In addition, these wall-function methods rely on care- 
ful integrations of the production and dissipation of the kinetic energy across the 
first mesh volume in order to have the mean kinetic energy identified with the mesh 
volume reflect the presence of the sublayer and buffer layers that exist there. In 
the first of the wall function methods (ref. 20), the "law of the wall" was used 
with fixed coefficients appropriate strictly for flat plates at constant pressure. 

In the second wall function method (ref. 21), these coefficients were allowed to 
vary, reflecting the behavior of the sublayer to strong axial pressure gradients, 
and in addition, provisions were made for the presence of the lower portion of the 
wake region of the boundary layer. These modifications are influential only in 
regions approaching separation, and in the separation zone, and immediately after 
reattachment. It should be noted that it was the added viscous terms (the so-called 
"low Reynolds number terms") which allowed integrations to the surface in the Jones 
Launder model that caused its difficulties and were eliminated through the use of 
the wall function approaches. 

Another turbulence model used in these comparisons was that of Wilcox and 
Rubesin (ref. 16). This model is similar to that of Jones and Launder in that two 
turbulence field equations also are employed to yield turbulence scale quantities on 
which to base the eddy viscosity. In this model, the variable used to develop the 
length scale is the specific rate of dissipation, i.e. the rate of dissipation 
divided by the kinetic energy. The turbulence model differs from the Jones-Launder 
model in another very significant way; it does not require additional viscous terms 
to permit integrations to the surface, but accomplishes this through the use of 
damped coefficients in the high Reynolds number terms. A drawback of this turbu- 
lence model is that the second equation uses the square of the specific dissipation 
rate as its dependent variable. This is a quantity that varies with distance from 
the surface to the fourth power, a variation that is too rapid for the second-order 
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finite differencing used in the computations to be accurate in the vicinity of the 
surface. To avoid inaccuracies, very near the surface an analytical expression is 
used to evaluate the dissipation rate variable. This analytical expression, how- 
ever, causes instabilities at the match point with the numerical solution of the 
second modeling equation. In addition, this model has particularly stiff mathemati- 
cal characteristics when wall functions are not used. These instabilities and 
stiffness characteristics prevented this model from being incorporated into the same 
algorithm as was used for all the other turbulence models (ref. 20). Consequently, 
the results presented here for the model of Wilcox and Rubesin were obtained with 
the older, less efficient, more complicated, but more robust hybrid code of 
MacCormack (ref. 17). The results of this model shown here are restricted to the 
surface pressure distribution, the only recoverable results available to the authors 
(ref. 15). 

To relieve the near-surface numerical difficulties associated with the Wilcox- 
Rubesin model, Coakley (refs. 6 and 22) developed turbulence models in which the 
dependent variables of the turbulence field equations were the slower varying square 
root of the kinetic energy and the first power of the specific dissipation rate. 
Coakley introduced low Reynolds number and damping terms that reduced or eliminated 
the types of unusual stiffness associated with the Wilcox-Rubesin and the Jones- 
Launder turbulence models and allowed efficient computation. In defining the coef- 
ficients of these new equations, Coakley transformed those of Jones and Launder into 
what is called here "Coakley-1" and those of Wilcox and Rubesin into "Coakley-2". 
Recall that the airfoil computations shown earlier were based on the Coakley- 1 
turbulence model . 


Results of Turbulence Model Comparison 

Figures 3(a) and 3(b) show the comparison of computed results based on the 
variety of turbulence models with the experimental surface-pressure data obtained in 
reference 11. The figures display the pressure coefficient distribution over the 
downstream half of the bump followed by one-half chord length of the cylinder. On 
these figures the experimental data, corresponding to a Mach number of 0.875 and a 
Reynolds number of 13.1 * xIO , are represented by the open circles, whereas the 
computed results are indicated with solid lines. The same data points are common to 
all the figures. Figure 3(a) compares the turbulence models of Cebeci-Smith, 

Johnson and King, and Coakley' s two versions of the two equation model. As 
expected, the Johnson-King model produces excellent results, in both locating the 
position of the shock wave and representing the pressure levels in the region of 
separation (see table 1). On the other hand, use of the classic Cebeci-Smith turbu- 
lence model yields a shock-wave location that is about 1/1 0th chord downstream from 
the experimental value, and much higher pressures in the region of separation. The 
results based on the use of the Coakley- 1 model are quite similar to those of the 
Cebeci-Smith model. This behavior illustrates why Coakley 's higher-order model 
results in the airfoil calculations shown in figure 1 were so similar to those from 
the Cebeci-Smith or Baldwin-Lomax algebraic models. The modified model, Coakley-2, 
shows some improvement in establishing both the location of the shock wave and in 
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the pressure level in the vicinity of separation. Figure 3(b) shows the same curves 
for the turbulence models of Jones-Launder , Wilcox-Rubesin , and three-examples of 
the Viegas-Rubesin wall function methods as applied to the Jones-Launder model. It 
is observed that the Jones-Launder model yields results quite similar to those of 
Coakley-2, with a slightly better prediction of the shock location but with essen- 
tially the same pressures in the separation zone. The other two-equation model, 
Wilcox-Rubesin, shows pressure coefficients that are uniformly lower than the 
others. These latter results are from computations made with the older hybrid code 
in which upstream boundary conditions were not consistent with requirements along 
characteristic lines as is made in the more modern codes. If the computed results 
from the Wilcox-Rubesin model are moved downward, i. e. the Cp are lowered uni- 
formly, they agree essentially with those from Jones-Launder or Coakley-2. These 
comparisons reveal that the two-equation turbulence models, except for Coakley-1, 
yield about the same pressure-distribution results regardless of dependent variable 
used in the second scale equation or of the differences in the details of how the 
equations are modified to permit integration to the surface. 

The remaining figures on figure 3(b) demonstrate the results from computations 
employing the wall-function techniques of references 20 and 21 with the Jones- 
Launder turbulence model. These techniques can reduce computation costs by about an 
order of magnitude (ref. 20). In these cases, the slightly low values of Cp shown 
ahead of the shock wave are believed to be caused by use of slightly high value of 
free-stream dynamic pressure. The basis of this observation is that comparisons 
directly in pressure show better agreement. 

The figure on the lower right represents the results from the simpler wall 
function version, reference 20. Two curves are presented, representing the results 
corresponding to different sizes of the initial mesh volumes adjacent to the sur- 
faces. These are identified by the value of y + , the dimensionless wall distance 
corresponding to the center of the mesh volume and values of the wall shear stress 
about a chord upstream of the bump's leading edge. Because an exponential version 
of the logarithmic "law of the wall" is used to define surface shear in this wall- 
function technique, it is desirable to assure that the y + lies well within the 
fully turbulent flow, that is, to have a value greater than about 30. Accordingly, 
two values of y + were tried, namely y + = 140 and y + = 60. It is observed that 
the computed location of the shock wave is sensitive to the choice of the y + . The 
pressure level in the separation zone, however, is not as sensitive to the 
different y + . The case with y + = 60 gives results that are a decided improvement 
over those from the original Jones-Launder model. For this wall-function scheme, it 
appears that the size of the first mesh cell adjacent to the surface becomes another 
modeling decision, or constant, in the process. Partially to remedy this matter, 
the second wall-function scheme was developed in reference 21. While reducing the 
sensitivity of the results to the value of y + , the effect was not entirely elimi- 
nated and the y + still remains a modeling parameter. The results with y + □ 140 
and the improved wall-function technique are shown in the lower-left of the fig- 
ure. These results are very good in establishing both the location of the shock 
wave and pressure levels in the region of separation, and are the nearest rival to 
the Johnson-King model. With y + = 60, for this wall-function model, values of 
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pressure coefficient indicated a movement of the predicted location of the shock 
wave of x/c = 0.02 from the y + = 140 case, and otherwise acted similarly to 
the y + = 60 case of wall function 1 . The y + = 60 case with the wall function 2 
model also showed some waviness in the pressure coefficient in the separated zone 
that we believe was caused by the absence of a buffer region in this model. 

When the values of skin friction evaluated by the various models are compared, 
it is found that they differ to a very large extent. Unfortunately, no measurements 
of skin friction were made in the experiment of reference 11 so that no conclusions 
regarding the relative merit of the the various turbulence models regarding skin 
friction can be made. The experiment did contain measurements of the mean points of 
separation and reattachment obtained by means of the oil film technique. These 
data and the points that bound the regions in the computations where skin friction 
is negative (i.e., the separation zone) are indicated in table 1. Again, the supe- 
riority of the Johnson-King model is demonstrated as it almost duplicates the exper- 
imental data. The wall-function techniques with chosen values of y + are next in 
performance, showing reattachment locations on the cylinder roughly 0.1 of the 
bump's chord upstream of the measured location. Finally, all of the two-equation 
models with integration to the surface show separation points approximately 0.3 
chord downstream of the measured values and reattachment points about 0.1 chord 
upstream. 

Figure 4 shows the development of the displacement thickness over the rear half 
of the bump and on the downstream cylinder. Recall that it is the displacement 
thickness that produces the influence of the viscous region upon the inviscid 
flow. As in figure 3, the experimental data are shown with open circles and are 
repeated on each of the graphs to facilitate the comparison of the various models. 

In figure 4(a), it is seen that the classic Cebeci-Smith model gives results that 
are quite inferior to the Johnson-King model. The Coakley-1 model, in this case, 
performs a bit better than does the Cebeci-Smith model, and the improvements in 
Coakley-2 are quite significant in the region of rising displacement thickness. In 
figure 4(b), the Jones-Launder model is shown again to yield results essentially the 
same as those of Coakley-2. The Viegas-Rubesin wall-function- 1 approach with 
y + = 140 produces similar results. The other wall-function approaches, with their 
selected values of y + , yield results for displacement thickness that surprisingly 
are even better than those of the Johnson-King model. 

Another interesting variable with which to compare the performance of the 
various turbulence models is the the shape factor. It is a quantity that enters the 
momentum integral equation used in some of the inviscid-viscous interactive schemes, 
and also can be used to infer the behavior of the momentum thickness, since all of 
the turbulence models yielded reasonably accurate values of the displacement thick- 
ness. Figure 5 shows the development of the shape factor over the latter half of 
the bump and on the downstream cylinder. One feature evident from all the computed 
results is a dip in the shape factor at the onset of separation, something that is 
not evident in the experimental data. The cause of this is unknown. Other than 
this, it is seen in figure 5(a) that the shape factor given by the Johnson-King 
model generally behaves the best of the models. The Cebeci-Smith and the Coakley-1 
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models fail to develop the steep rise in the shape factor evident in the data in the 
separation region. This is evidence that the momentum integral computed by these 
models in this region must be much higher than the actual. This behavior is cor- 
rected to some extent by the Coakley-2 model and the Jones-Launder model shown in 
figure 5(b) The wall function models give mixed results. The case of the Viegas- 
Rubesin wall function 1 with y + = 60, which gave the better results for the other 
variables considered previously than did the case with y + = 1 40 , overshoots the 
data in the region just aft of the point of separation. This shows that the local 
momentum thickness computed with y + = 60 must be yielding low values of the momen- 
tum thickness there. The Viegas-Rubesin wall-function-2 computation with y + = 140 
shows similar tendencies but behaves quite well over the downstream cylinder. The 
irregularities shown in this latter figure are believed to be the result of locally 
low values of the y + in the region of separation. This can allow the center of 
the first mesh volume off the surface to be either within a fully turbulent and/or 
fully laminar zone of the boundary layer (one of the features that is permitted by 
this version of the wall function models). Blending these distinct layers with a 
buffer layer could possibly have avoided the bumps and valleys of the computation. 

To further illustrate the relative behavior of the various turbulence models, 
the velocities in the axial direction predicted by the computations are compared in 
figure 6 with the corresponding measured data. All the profiles are at stations 
within the experimentally determined separation zone over the trailing-edge region 
of the bump and the cylinder downstream. The velocities are normalized by the free- 
stream velocity and plotted as functions of the vertical distance from the surface 
expressed in centimeters. With linear coordinates, a plot such as this only shows 
the behavior of what is usually termed the outer portion of the boundary layer. 

The sublayer and buffer layer cannot be seen distinctly without recourse to loga- 
rithmic plotting. These profiles are shown at values of x/c = 0.75, 0.875, 1.00 
(the trailing edge of the bump), and 1.125. Recall that separation and reattachment 
was indicated at x/c = 0.67 and 1.17. At x/c = 0.75, roughly 1/10th chord down- 
stream of the mean separation, it is surprising that the Johnson-King model does not 
show the best agreement with the velocity data in view of its success in predicting 
the surface pressure and integral thicknesses. A possible explanation of this 
behavior is given in reference 12 where it is shown that better agreement between 
the velocity computations and the experimental data can be achieved by shifting the 
computed results downstream by only about 0.03 chord. Another way of arguing this 
point is to note that at x/c = 0.75 very rapid changes to the flow field are 
occurring so that small misalignments between the computations and the experimental 
data in the streamwise direction cause large apparent profile errors. Subject to 
this observation as a caveat, it is found that the Coakley-2 and Jones-Launder 
models give the best fits to the velocity data at this station. Since the flows for 
these latter models have not yet separated at this station, their agreement with the 
far field velocity data shows how sensitive near-wall turbulence modeling is to the 
prediction of surface phenomena. Again, there are similarities in the results of 
the Cebeci-Smith and Coakley-1 models, and in the wall-function models 1 and 2 
with y + = 60 and 140, respectively. Errors in the velocities toward the outer 
part of the boundary layer are indications of poor representation of the surface 
pressure at a particular station, and indicate the largest errors for the 
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Cebeci-Smith and Coakley-1 models. This observation is consistent with the pressure 
coefficient results shown on figure 3(a) at x/c = 0.75. 

The remaining velocity profiles are at stations where flow-field changes are 
less rapid, and are more indicative of the performance of the various turbulence 
models. Except for the Cebeci-Smith and Coakley-1 models, all the turbulence models 
give results for the velocity profiles that are very representative of the experi- 
mental data. Apparently, far-field velocity data in and near a region of separation 
are not critically sensitive to turbulence models employed. Of course, in the near- 
wall region, this is not the case as the skin friction demonstrated by the various 
models differed considerably in magnitude and even in sign at some stations. Also, 
the position of separation is sensitive to the near-wall velocity behavior, and also 
differed considerably with the different turbulence models. 

A final comparison of the performance of the various models is shown in fig- 
ure 7, where profiles of the Reynolds shear stress parallel to the streamwise direc- 
tion are compared with the corresponding experimental measurements. Recall that a 
linear plot in distance from the surface shows principally the behavior of the outer 
part of the boundary layer. The stations shown cover the experimental region of 
separation and are the same as those used for the mean velocities in figure 6. 

The most striking feature of the series of plots on figure 7 is the similarity 
of the general behavior all the turbulence models. This behavior is rather surpris- 
ing when consideration is given to fact that the zones of separation predicted by 
many of the models had little relationship to the experimental zone covered by these 
figures. For example, the Coakley-1 model barely shows any separation, yet the 
Reynolds stresses predicted by it are not that grossly different from those of the 
Johnson-King model. From this it can be inferred that the Reynolds stresses outside 
the viscous sublayer and buffer region are rather insensitive to the presence or 
absence of separation. This observation is consistent with the concept in the 
algebraic eddy viscosity models that the the proper scaling in the outer part of the 
boundary layer is the displacement thickness and the velocity at the boundary layer 
edge, which are more similar in the computations based on the various models than 
are their zones of separation. 

In detail, it is observed that none of the models yields values of Reynolds 
shear stress as large at the measured values at the downstream stations. The 
Cebeci-Smith and Coakley-1 models perform the poorest in this regard. The more 
successful models, however, are those whose peaks of shear stress are larger in 
magnitude and occur closer to the distances from the surface shown by the maxima in 
the experimental data. Modeling the outer region of the boundary layer seems to be 
important in this flow, and may explain why models that account for the out of 
equilibrium conditions there, through a lag equation or field equations, are the 
more successful. 
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CONCLUDING REMARKS 


It is generally accepted that there is no unique way of modeling turbulence in 
transonic flow fields. Consequently, a model can only be gauged by its degree of 
success in computing a particular type of flow field. To do this, it is necessary 
to eliminate, as much as it is possible to do so, experimental errors in the data 
used for model assessment. We tried to do this here by dealing with a two- 
dimensional flow field that was free of uncertainties due to wind tunnel wall 
effects and to unknown angles of attack, in contrast to the ambiguities that exist 
in the current crop of airfoil experiments that were used to assess computation 
codes and turbulence models earlier. Also, past assessments dealt only with a 
limited number of turbulence models; two algebraic eddy viscosity models, a lag 
model that determines the growth of the maximum shear stress along a surface that is 
then used to scale an algebraic eddy viscosity across the boundary layer, and a 
single version of two-equation modeling where the variables upon which an eddy 
viscosity is evaluated are found from field equations. In the present paper, five 
more two-equation eddy viscosity models have been considered, two of which involve 
the use of wall functions to represent surface boundary conditions and avoid carry- 
ing integration to the surface. 

The preceding descriptions of the comparative performance of these turbulence 
models reflect the two main underlying turbulence modeling philosophies. One 
approach argues that for a specific class of flows it is best to begin with a simple 
model and adapt it just enough to capture the most essential elements of the new 
flow. The other approach is to apply a so-called "universal" model. These employ 
field equations for their basic variables and can be applied without change to any 
type of flow. In principle, a good version of the latter type of model would yield 
accurate quantitative results for a disparate variety of flows. In practice, how- 
ever, it is found that while the models can be applied rather easily computationally 
to new situations, the results they give in particular cases are often not as accu- 
rate as would be desired for aerodynamic design purposes. 

The Johnson-King model is representative of the first approach applied to the 
development of an accurate model for the prediction of transonic flows over simple 
aerodynamic bodies. A key element in its development was the recognition that it 
was essential to account for the nonequilibrium between the turbulence and mean flow 
in the outer portions of the viscous layers. In contrast, the earlier models of 
Cebeci-Smith and Baldwin-Lomax assume complete equilibrium between the turbulence 
and the mean motion. Johnson and King introduced nonequilibrium in the streamwise 
direction by utilizing an equation to evaluate the rate of development of the maxi- 
mum shear stress in the boundary layer. This maximum shear stress was then used to 
scale the eddy viscosity in algebraic expressions across the boundary layer. This 
latter step implies an equilibrium in the direction normal to the surface. The 
model required the introduction of two additional modeling constants beyond those 
used in the algebraic eddy viscosity models. As was shown earlier, this Johnson- 
King model could be adjusted to produce excellent results for the mean velocities 
and surface pressures over the axisymmetric bump. This is not surprising as these 


595 



experimental data were used to define the additional modeling constants used in the 
model. Significantly, however, it was found that this model gave good comparisons 
with low-speed data (ref. 23) and also, as was demonstrated earlier, generally 
produced the best comparisons among the Navier-Stokes methods with the three sets of 
independent airfoil data. This suggests that the Johnson-King model has the breadth 
to be applied with confidence in a predictive sense to some new transonic flow 
provided it has separation zones comparable in size to those in the flow over the 
transonic bump. 

Several examples were shown of the performance of "universal" two-equation eddy 
viscosity models integrable to the surface (refs. 6, 19, and 22). These models 
utilized coefficients that were based on flows in pipes or attached boundary layers, 
and were unmodified to reflect the requirements of transonic flows with shock 
waves. Although these models allow the turbulence and mean motion to be out of 
equilibrium in both the streamwise and normal directions, it was found that the 
Coakley-1 model gave essentially the same overall results as the equilibrium model 
of Cebeci-Smith. Note that computations of the transonic flow over the bump have 
not been made with the Baldwin-Lomax model, but, based on past modeling experience 
(ref. 2) this model is expected to yield results similar to those of the Cebeci- 
Smith model. Some improvement occurred with the Coakley-2 model, where the modeling 
coefficients were altered from those of Coakley-1, but the results were no better 
than given by the earlier Jones-Launder model. It is not clear why the Coakley-2 
model, designed to account more accurately for separated flows at supersonic speeds, 
did not give more accurate results for the transonic flow considered here. Perhaps 
the level of effort of trying to improve the Coakley-1 model was less than was 
exerted by Johnson and King; or the act of improving the field equation models is 
just inherently more difficult than with the simpler models. In either case, the 
Coakley-2 model still leaves something to be desired in the prediction of transonic 
flows . 

The other models shown were applications of the use of wall functions to the 
Jones-Launder model. Implicit in the use of wall functions is the assumption that 
the near-wall region is one where equilibrium between turbulence and mean flow 
exists. These wall functions were developed principally to save computer effort by 
reducing the numerical stiffness and the number of mesh points employed. They have 
not been tailored for any particular flow or two-equation turbulence model and are 
quite general. The wall functions have been adapted to three two-equation models 
of turbulence (two k, epsilon models and the Wilcox-Rubesin model) and successfully 
applied to a variety of complex flows and over a range of Mach numbers. Generally, 
the wall functions have made the computer codes more robust and, with a judicious 
choice of the size of the first mesh volume, can obtain enhanced agreement between 
the numerical results and the experimental measurements. Such improvements had 
previously been observed in several transonic duct flow calulations (refs. 20 
and 21). Their use here, with dimensions of the first mesh volume that resulted in 
values of y + = 60 and 140 for Viegas-Rubesin wall functions methods 1 and 2, 
respectively, also reduced computation costs and improved the comparison of the 
computations with the experimental results of the transonic bump flow. It should be 
noted that these wall functions were of the type that perform accurate integrations 
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of the production and dissipation of the turbulence kinetic energy within the first 
mesh volume off the surface. This allows accounting for the larger fraction of the 
viscous-flow region compared to fully turbulent flow region within this mesh volume 
that occurs at the stations in or near the separated flow. 

It can be concluded from this review that algebraic eddy viscosity turbulence 
models that assume complete equilibrium between the turbulence and the mean motion, 
when coupled with the Navier-Stokes equations, are inadequate to yield the accura- 
cies in the calculation of transonic flows that are necessary for aerodynamic 
design. Based on comparisons with airfoil data and over an axisymmetric bump, it 
has been found that the Johnson-King model remedied this situation in two dimen- 
sional flows by merely adding the complexity of an additional ordinary differential 
equation to account for the non-equilibrium growth of the maximum shear stress in 
the boundary layer along the surface. Equilibrium of the turbulence and mean flow 
is still assumed across the boundary layer. Of course, if ap airfoil design con- 
tains a surface jet, or multiple control surfaces, situations can develop within the 
boundary layer where more than one maximum in shear stress can occur and ambiguity 
in the model develops. Also, the extension of the concept of basing algebraic eddy 
viscosities on a maximum shear stress to three-dimensional flows will require con- 
siderable development and loss of simplicity. 

Two-equation models that can be integrated to the surface require further 
careful development before they can be considered to be predictive models. Emphasis 
in this development should be on the improvement of performance in the near-wall 
regions under severe adverse streamwise pressure gradients and in regions of separa- 
tion. The effort to do this is worthwhile because of the particular adaptability of 
this class of closure models to new situations. Some examples of this are: three- 
dimensionality, surface jets, and multi-element airfoils and wings. In design 
computations involving many changes of parameters, the cost savings introduced by 
wall functions should not be overlooked. Thus, wall functions should be improved to 
further reduce the sensitivity of the results to differences in the dimensions of 
the first mesh volume off the surface. 

It is also significant that, to date, some of innovations to turbulence model- 
ing that have shown some promise in low-speed, incompressible flows have not been 
applied yet to transonic flows. One example of this is the algebraic stress model 
that eliminates the controversial eddy viscosity constitutive relationship 
(ref. 24). This approach has inherent advantages when isotropies in the Reynolds 
stresses are expected to depart from those associated with two-dimensional shear 
flows. Examples of situations where this is likely to happen are three-dimensional 
flows, regions of intense streamwise curvature, and in the presence of shock 
waves. Another example of developments in incompressible flow computations is the 
four-equation approach, where the shear stress, the most important Reynolds stress, 
is determined through an additional field equation and the different roles played by 
the large and small eddies of turbulence, which primarily exchange momentum and 
dissipate turbulence, respectively are introduced through another kinetic energy 
equation (ref. 25). In addition, it has been shown that modeling which represents 
each Reynolds stress with a field equation and utilizes a single-scale equation can 
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be applied at acceptable computer costs to transonic flows (ref. 26), and fine- 
tuning of this approach to improve its accuracy is warranted not necessarily as a 
design tool, but as a standard of comparison and guide for the simpler approaches. 
Finally, in view of the success of some of the interactive inviscid/viscous methods 
in predicting the airfoil lift and drag coefficients of figure 1, it would be inter- 
esting to see these methods applied to the transonic bump flow and to compare the 
results to those presented here that were based on solutions to the Reynolds- 
averaged Navier-Stokes equations. 
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TABLE 1.- LOCATION OF SEPARATION AND REATTACHMENT, x/c 



Separation Point 

Reattachment Point 

Experiment 

0.67 

1.17 

Cebeci-Smith 

0.73 

1.04 

Johnson-King 

0.68 

1.17 

Coakley 1 

0.99 

1 .00 

Coakley 2 

0.96 

1 .04 

Jones-Launder 

0.94 

1.04 

Viegas-Rubesin 
Wall Function 1 
y + = 60.0 

0.68 

1 .08 

y + = 140.0 

0.85 

1.06 

Viegas-Rubesin 
Wall Function 2 
y + = 140.0 

0.69 

1.11 


I 


l 
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Figure 1.- Comparison of measured and computed force coefficients of transonic 
airfoils. Symbol legend: experimental data •; solutions from interactive 

methods O; Navier-Stokes solutions with Baldwin-Lomax model □; Cebeci-Smith 
model Pj Johnson-King model <£>; and Coakley-1 model A- 
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(c) VA-2, Re c = 6 x 10 6 . 
Figure 1.- Concluded. 


Re 00 = 13.1 X 10 6 /m 



Figure 2.- Schematic diagram of transonic bump model and wind tunnel conditions 

(ref. 11). 
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(a) Turbulence models of Cebeci-Smith, Johnson-King, Coakley-1, and Coakley-2. 



x/c x/c 

(b) Turbulence models of Jones-Launder , Wilcox-Rubesin , Viegas-Rubesin 
Wall Functions 1 and 2. 


Figure 3.- Pressure distribution over the circular arc bump. The same experimental 
data, represented by open circles, are repeated in each figure. 
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(a) Turbulence models of Cebeci-Smith, Johnson-King, Coakley-1, and Coakley-2. 



(b) Turbulence models of Jones-Launder , Wilcox-Rubesin, Viegas-Rubesin 
Wall Functions 1 and 2. 

Figure 5.- Shape factor over circular arc bump. The same experimental data, 
represented by open circles, are repeated in each figure. 
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(a) Turbulence models of Cebeci-Smith, Johnson-King , Coakley-1, and Coakley-2. 


Figure 6.- Mean velocity profiles over the separation zone on the circular arc 
bump. The same experimental data, represented by open circles, are repeated in 
each figure. 
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(b) Turbulence models of Jones-Launder , Wilcox-Rubesin, Viegas-Rubesin 
Wall Functions 1 and 2. 


Figure 6.- Concluded. 



I 


X/C = 0.75 X/C = 0.875 X/C = 1.000 X/C = 1 .125 



(a) Turbulence models of Cebeci-Smith, Johnson-King, Coakley-1, and Coakley-2. 

Figure 7.- Shear stress profiles in the separation zone on the circular arc bump. 
The same experimental data, represented by open circles, are repeated in each 
figure. 
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X/C = 0.75 



X/C = 0.875 X/C = 1.000 X/C = 1.125 







-U'V'/U 2 ^ -U'V'/U 2 *, -U'V'/U 2 *, 


(b) 


(b) Turbulence models of Jones-Launder , Wilcox-Rubesin, Viegas-Rubesin 
Wall Functions 1 and 2. 


Figure 7.- Concluded. 
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